On the Nagaoka polaron in the t — J model 
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It is widely believed that a single hole in the two (or three) dimensional t — J model, for sufficiently 
small exchange coupling J, creates a ferromagnetic bubble around itself, a finite J remnant of the 
ferromagnetic groundstate at J = (the infinite U Hubbard model), first established by Nagaoka 
0. We investigate this phenomenon in two dimensions using the density matrix renormalization 
group, for system sizes up to 9 x 9. We find that the polaron forms for J/t < 0.02 — 0.03 (a somewhat 
larger value than estimated previously). Although finite-size effects appear large, our data seems 
consistent with the expected 1.1( J/i)~ 1//4 variation of polarion radius. We also test the Brinkman- 
Rice model of non-retracing paths in a Neel background, showing that it is quite accurate, at larger 
J. Results are also presented in the case where the Heisenberg interaction is dropped (the t — J z 
model). Finally we discuss a "dressed polaron" picture in which the hole propagates freely inside a 
finite region but makes only self-retracing excursions outside this region. 



An understanding of the groundstate of the two dimen- 
sional t— J model near half-filling remains one of the most 
urgent open problems in condensed matter theory. The 
Hamiltonian may be written 



H = 



E 

<i,j> 



-tP(J2 4a c j<* + h.c.)P + J{S, ■ S - riirij/4) 



(1) 

Here h.c. stands for Hermitean conjugate, Ci a annihilates 
an electron of spin a at site i, P projects out states with 
no double occupancy on any site, < i,j > represents 
pairs of nearest neighbor sites, and: 



Si 



0*0/3 
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Even the case of a single hole relative to half-filling is 
highly non-trivial. At J = 0, corresponding to the infi- 
nite U Hubbard model, the groundstate was proved to 
be ferromagnetic by Nagaoka. jj] This happens because 
the hole can move around most efficiently, thus minimiz- 
ing the kinetic energy, when the spins are fully polarized. 
Conversely, an antiferromagnetic exchange coupling, J, 
favours anti-parallel spins leading to the picture that a 
bubble of polarized spins forms around the hole while the 
spins further away are in the antiferromagnetic ground- 
state. The hole essentially moves around freely only in- 
side this bubble, or polaron, longer range motion leading 
to disruption of the antiferromagnetic groundstate. How- 
ever, this picture is only expected to hold, at best, for 
sufficiently small J/t. At somewhat larger J/t quite dif- 
ferent behavior is expected, perhaps the Brinkman-Ricc 
picture in which the hole makes only self-retracing ex- 
cursions from the origin in order to avoid disrupting the 
antiferromagnetic background. At values of J relevant 
to the cuprates, J/t ~ 0.3 — 0.4, numerical studies have 
demonstrated the presence of frustrating singlet correla- 
tions across mobile holes, connecting sites on the same 



sublattice |2j. However, a simple analytic approach re- 
producing this effect is lacking. 

Recent advances in numerical methods now make it 
possible to study the one hole problem accurately. We 
present here density matrix renormalization group Q| 
(DMRG) results for the small J/t regime on lattices of 
' size up to 9 x 9 with boundary conditions that are ei- 
ther open or mixed open and periodic (cylindrical). The 
number of states kept per block ranged from 2000 to 
4000. Various methods of studying the Nagaoka polaron 
are possible. We study the nearest neighbour exchange 
energy (Sp ■ SV+e — npnp + ^/A), as a function of position 
r. This is increased in the vicinity of the hole, which 
resides near the center of the lattice due to a combi- 
nation of boundary conditions and initial conditions in 
the DMRG calculation. We calculate the groundstate 
energy as a function of J/t, showing reasonable agree- 
ment with the y/j dependence expected in the polaronic 
state, reviewed below. Most directly, we can calculate 
the total spin of the groundstate with an even number of 
sites and a single hole. We find that this varies from 1/2 
for J/t > 0.02 — 0.03, where the polaron disappears, to 
(N— l)/2, where N is the number of sites for J/t ~ .001. 
In this latter case the polaron becomes larger than the 
lattice size. 

The energy and radius of the Nagaoka polaron, as a 
function of J/t are derived by assuming that the hole 
wave-function is that of a free particle with vanishing 
boundary conditions at the edge of the circular bubble. 
The groundstate is given by the lowest energy such state. 
R is determined by balancing the kinetic energy of the 
hole against the energy cost of the ferromagnetic bonds. 
For small J and hence large R the hole wave-function 
corresponds to a small wave- vector, at the bottom of the 
band, which can be approximated as quadratic, with an 
effective mass of l/2t. (Henceforth we set t=l.) The 
wave- function is hence a Bessel function, Jo(kr) with k 
chosen so that Jo(kR) — 0. The first zero of Jq(x) occurs 
at x Rj 2.40. Hence the kinetic energy, compared to that 
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at half-filling is: 

E K = -4 + {2A0f/R 2 . (3) 

To this we must add the energy cost of the ferromagnetic 
bonds. A ferromagnetic bond costs energy J/2 relative to 
the energy of a bond in an up-down state. However, the 
actual energy per bond of the Neel groundstate for the 
two dimensional square lattice is approximately — 0.33J, 
so the energy cost of a ferromagnetic bond is about 0.58 J. 
Taking into account that there are 2 bonds per unit area 
(we set the lattice constant to 1) we see that the magnetic 
energy, measured relative to that of the groundstate with 
no hole, is: 

E M = l.l6JnR 2 . (4) 

Thus the total energy of the polaron, relative to the an- 
tiferromagnetic groundstate with no hole present is: 

E = -4+ {2AQ) 2 /R 2 + 1.16J7ri? 2 . (5) 

Minimizing the energy with respect to R gives: 

R = 1.12J~ 1/4 

E = -4 + 9.20\/j. (6) 

The total spin of the polaron is: 

S = (1/2)ttR 2 = 1.97 J- 1/2 . (7) 

We note that, even if we assume that these formulas 
become exact in the limit J — ► for an infinite area sys- 
tem, we should expect severe finite size effects. R must 
be much bigger than 1 lattice spacing in order for the 
assumption of a circular polaron with a simple vanish- 
ing boundary condition to be reasonable. On the other 
hand, the system length must be much larger than R in 
order for the polaron not to be perturbed by the finite 
system size. Thus, in practice there will be at best a 
narrow range of small J's for which these formulas apply 
in a numerical simulation. If J is too small, the polaron 
simply fills the entire system. In this case its energy is 
trivial to calculate. The kinetic energy is that of a sin- 
gle hole in a filled (spin polarized) band. If the system 
has dimensions L x M with free boundary conditions, for 
example, this is: 

E K = -2cos[7r/(£ + l)] - 2cos[tt/(M+ 1)]. (8) 

(Si ■ Sj — riirij/A) is exactly zero in this state. The mag- 
netic energy, relative to the no-hole ground state, is sim- 
ply 0.58 times the number of bonds 

E M = 0.58 J[2LM - (L + M% (9) 

up to finite size corrections to the groundstate energy of 
the hole-free system. 



We expect the polaron to fill the entire system when 
J is reduced to a value of order l/N 2 where N is the 
number of sites. On the other hand, as J is increased 
we expect the polaron to begin showing departures from 
circular shape as R becomes not so much larger than 1. 
Eventually when R is of 0(1) the magnetization of the 
polaron is reduced to 1/2 at which point it no longer ex- 
ists. (S = 1/2 is the smallest possible value for one hole 
relative to half-filling in a system with an even number 
of sites.) Thus, as we decrease J for a given large system 
size we expect S to increase approximately as 1.97J -1 / 2 
eventually saturating at (N — l)/2 while the energy de- 
creases as —4 + 9.2 J 1 / 2 until it eventually crosses over 
to the linear J-dependence of Eqs. (||) and (||). An- 
other way that the polaron can be studied using DMRG 
takes advantage of the iterative nature of this numerical 
method. We begin our DMRG sweeps with the hole sit- 
ting at the origin. In principle, after an infinite number 
of DMRG sweeps in a periodic system, it should go into 
an extended state where it spreads over the entire lat- 
tice. In this case the groundstate should be a momentum 
eigenstate in which observables like the magnetic energy 
density are translationally invariant. However, in prac- 
tice, we use open or cylindrical boundary conditions. For 
open BCs, the boundaries tend to localize the polaron in 
the center of the system. Even with cylindrical BCs, we 
find that the polaron tends to remain localized at the 
origin. The fact that this occurs, while the expectation 
value of the Hamiltonian in the approximate groundstate 
appears to be close to its minimum, indicates that the 
polaron has a large effective mass. Since the polaron re- 
mains well localized we can study its size by measuring 
spin correlations, (S?- Sf+e function of r. 

These measures of the polaron size versus J are shown 
in the figures. In Fig. 1 we show the kinetic and mag- 
netic energy versus J. Note that the agreement with Eqs. 
d) and (§J) is fair for J < 0.02 - 0.03 but that the data 
starts to deviate strongly from the theoretical curve at 
this point, suggesting the breakdown of the polaron pic- 
ture. The sharp features in Em and Ek at J rs 0.03 are 
rather surprising and may be numerical artifacts. The 
uncertainties in measurements of the magnetic and ki- 
netic energies separately are much larger than those as- 
sociated with the total energy. We have included rough 
error bars, estimated using the energy measurements as 
a function of the number of states kept per block, and 
have also shown some data for a 9 x 9 lattice for this 
reason. The 9x9 results do not clearly show a local min- 
imum in Em near 0.03. Regardless of whether there is 
such nonmonotonic behavior in Em or Ek, there seems 
to be rather sharp behavior associated with an abrupt 
crossover as the polaron disappears. We have also mea- 
sured the spin of the polaron for various J, shown in Fig. 
2. While S certainly increases with decreasing J, even- 
tually reaching saturation, at J w 0.001, it seems to lie 
well below the polaron prediction of Eq. (pi) , at all J. 



2 



Finally, in Fig. 3 we present spin correlations, showing 
a disturbance of the antiferromagnetic state in the vicin- 
ity of the hole, situated at the origin. The region over 
which a significant departure from the antiferromagnetic 
state occurs gives another measure of the size of the po- 
laron. The polaron does indeed look approximately cir- 
cular, as expected. 

Taken together, we think that this data suggests rather 
strongly that the Nagaoka polaron does indeed form, and 
that Eqs. (^) and (|4|) are asymptotically correct, for small 
enough J and an infinite lattice. We have also shown that 
the value of J where this polaron picture breaks down 
is about 0.02 — 0.03, five times larger than a previous 
analytical estimate ||. However, this value is still about 
ten times smaller than that believed to be relevant to the 
cuprates. 
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FIG. 1. (a) The magnetic energy Em of a single hole as a 
function of J. The continuous line is Eq. (4), using Eq. (6). 
(b) The kinetic energy Ek- Here, the continuous line is Eq. 
(3), using Eq. (6). 
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FIG. 2. The total spin of the ground state of an 8 x 8 
system with one hole, as a function of J. The spin was deter- 
mined by looking for degeneracy in the ground state energy 
for different DMRG calculations as we very S z - The result 
for S is the largest value of S z degenerate with S z = 1/2. 
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FIG. 3. The expectation value of Si ■ Sj — mnj / 4, displayed 
as the thickness of each line, for a 9 x 9 open system with 
one hole and J = 0.01. The ferromagnetic polaron is clearly 
visible in the central region, where the expecation value is 
nearly zero. 

We have also compared the energy to that of the 
Brinkman-Rice approximation, [Q] in the form used by 
Shraiman and Siggia [|| (BRSS approximation). This 
approximation should be better for the t — J z model, in 
which the spin-flip part of the Heisenberg interaction is 
dropped. Then as the hole propagates it leaves behind a 
string of flipped spins. If the hole only propagates along 
the x axis, for example, then after it has moved r-steps 
there are (2r + 1) ferromagnetic bonds so the exchange 
energy is simply: 



V(r) = (5/2 + r) J, (r > 1) 
V(0) = 2 J. 



(10) 
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(Note that we include the —niUj/A term in the exchange 
energy.) In the BRSS approximation, the potential en- 
ergy is assumed to have this value for all strings of length 
r, regardless of their shape. However, if the hole moves 
around a single plaquette, then the exchange energy is 
only 5J after three steps, rather than (5/2 + 3) J. Thus 
the potential energy is not the same for all strings of 
length 3 and larger. Another aspect of the BRSS ap- 
proximation is to treat each string as a distinct state. 
In this case the problem becomes equivalent to a single 
partice on a Bethe lattice with a hopping term —1 and 
the potential energy of Eq. ([l0|), for all sites at depth r. 
The hole can hop in any of 4 directions at the first step. 
Thereafter, there are only three non-retracing directions 
available. Therefore, the Bethe lattice has 1 state at 
depth 0, 4 states at depth 1 and 4 x 3 r_1 states at depth 
r for r > 2. However, the Bethe lattice approximation is 
not strictly correct, since there are several different paths 
which lead to the identical hole position and spin orien- 
tation. For instance, the same configuration is obtained 
when the hole goes around a single plaquette 1-1/2 times 
either in a clockwise or counter-clockwise sense. Thus 
the Bethe lattice doesn't have the correct connectivity. 
The actual t — J z model would correspond to a Bethe 
lattice with some of the sites identified, starting at depth 
6, and potential energies which are not the same for all 
sites at a given depth. 

It is generally assumed, in the BRSS approximation, 
that the hole wave- function has the same value at all sites 
of a given depth. Then the Bethe lattice model becomes 
equivalent to a one-dimensional model with hopping am- 
plitude —2 between site and 1 and — \/3 between sites 
r and r+1 for r > 1. The potential energy is 2 J on site 
and (5/2+r) J on site r for r > 1. A free boundary condi- 
tion occurs at r = 0. The single particle groundstate for 
this one dimensional problem can easily be found numeri- 
cally. In particular, when J <C t, we can use a continuum 
approximation, approximating the dispersion relation as 
quadratic near the bottom of the band, giving: 

d 2 

H = -2V3 - V3— T + Jx, (a? > 0) (11) 
ax 

with a vanishing boundary condition at x = 0. An exact 
scaling argument then implies that the groundstate en- 
ergy is proportional to J 2 / 3 and the numerical solution 
of this eigenvalue problem (in terms of the Airy function) 
gives the energy: 

E « -2V3 + 2.81J 2/3 . (12) 

For larger J the continuum approximation breaks down 
but the model can be solved numerically on the lattice. 
Shraiman and Siggia H state that the numerical solution 
of this equation is well approximated by Eq. ( |l2| ) with 
the factor in the second term replaced by 2.74 over the 
range 0.005 < J/t < 1. (We have verified this result.) In 



Fig. 4 we compare the DMRG results for the groundstate 
energy to the BRSS approximation. The agreement looks 
quite good for J > 0.02 - 0.03. We note that the BRSS 
approximation does not give a lower bound on the energy 
and appears to lie below the actual energy as determined 
by DMRG. 




J 

FIG. 4. The total energy of a hole, measured relative to 
the undoped system, as a function of J, on a 10 x 8 system. 
The curve labeled BRSS is Eq. (12) with 2.81 replaced by 
2.74, and with an extra term J added on to account for the 
term —JniTij/4: in the Hamiltonian which was not included in 
the Hamiltonian studied by BRSS. The curve labeled Nagaoka 
is Eq. (6), also with an extra J added on. 

Some insight into the corrections to the polaron pic- 
ture for slightly larger J can be obtained by combining 
it with the BRSS string approximation. We could thus 
consider a state where the hole moves in a ferromagnetic 
background inside a bubble of radius R but also makes 
excursions outside this bubble which take the form of 
self-retracing walks or strings. Following BRSS we may 
approximate the lattice outside the polaron as a Bethe 
lattice. In other words, as shown in Fig. 5, we effectively 
attach "Bethe hair" to each point on the boundary of the 
polaron so that the hole propagates on a square lattice 
inside the polaron but on a Bethe lattice outside the po- 
laron. For large R (small J) the energy is close to —4. 
This implies than when the hole enters the Bethe lattice 
it propagates below the band-edge at — 2V3. Again mak- 
ing the reasonable approximation that the wave-function 
has the same value at all points at the same depth on 
the Bethe lattice, the problem is again approximated as 
a one-dimensional tight-binding model with a linear po- 
tential. The linear potential plays little role at small J 
since the sub-band wave-function decays exponentially 
even at J = 0. To see this note that the Schrodinger 
equation, everywhere on the Bethe lattice except near 
the top has the form: 

E0 r = ~V3(0 r + r+ i), (13) 
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with E w —4. The two solutions of this equation 
asymptotically behave as 4> r ~ 3 ±r / 2 . A valid solution 
of the eigenvalue problem must give only the exponen- 
tially decreasing solution, oc 3~ r / 2 . Actually solving 
the Schroedinger equation would require making small 
independent adjustments in the energy (near —4) and 
adjusting the wave-function near the edge of the polaron 
to ensure that only the exponentially decreasing solu- 
tion is obtained on the Bethe lattice. We might expect 
the Bethe lattice approximation to actually work better 
in this regime (smaller J when the energy approaches 
—4) than at larger J where the BRSS approximaton is 
normally applied (E > — 2v3), because the hole has 
only a small amplitude to descend to large depths in 
the Bethe lattice, due to the exponential decay of the 
wave- function. As noted above, the Bethe lattice approx- 
imation only starts to break down at fairly large depth 
(depth 6 for the lattice connectivity). We note that the 
Bethe hair may be regarded as a sort of boundary layer 
attached to the polaron. By the above argument it only 
has a width of a few lattice spacings, independent of R, 
at large R. Thus the naive estimates of R(J) and E(J) in 
Eq. (^|) should remain true at small enough J. The Bethe 
hair just represents a sort of 1/R correction to the po- 
laron approximation. However, we see that it will make 
the finite size difficulties with numerical simulations even 
more severe. The lattice must be large enough to acco- 
modate not only the polaron but also the boundary layer. 





spins inside the polaron are oriented in the xy plane, 
rather than along the z-axis. Classically this costs an 
exchange energy of J/4 per bond, rather than J/2 which 
would occur for z-orientation 0. Eqs. (||), @ are thus 
modified to: 



E M = (l/2)J z irR 2 
R = -1.39J7 1/4 



E = -4 + 6.03 V J 



(14) 



Our numerical results are consistent with this picture; we 
find S z = for all J. A comparison of the polaron and 
BRSS energies to our DMRG results in shown in Fig. 6. 
We also show the formula of Barnes, et al., 



E 



-3.63 + 2.93 J 



2/3 



(15) 



obtained from fitting Lanczos and Monte Carlo data. 
Perhaps surprisingly, the agreement with BRSS at larger 
J is a bit worse than in the t — J model. We also present 
(SpSZ +i — npnp +& /4:) in Fig. 7. A central region where 
this is close to —1/4 is visible, corresponding to a po- 
laron with spins lying in the xy plane. The slightly worse 
agreement with the polaron picture in the t — J z model 
may just reflect the fact that the boundary layer is wider 
leading to larger finite size corrections. A recent exten- 
sive discussion of the t — J z model may be found in Ref. 




FIG. 5. A schematic diagram of a polaron with Bethe 
lattice "hair". 

We also present some results on the t — J z model. 
One might again expect the polaron picture to hold at 
small enough J z . However, we no longer have such a 
direct measurement of the polaron size as was afforded 
in the Heisenberg case by the total spin of the ground- 
state, since S is not conserved in the Ising case. We can 
still measure the total z-component, S z . However, it is 
reasonable to expect that the ferromagnetically aligned 
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FIG. 6. For the t — J z model, the total energy of a hole, 
measured relative to the undoped system, as a function of 
J, on a 9 x 9 system. For comparison, we also show the 
approximate energies for a Nagaoka ferromagnetic polaron 
and for a hole according to the BRSS treatment. The same 
extra term J discussed in the caption of Fig. 4 was added to 
each of the three analytic curves. 
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FIG. 7. For the t — J z model, the expectation value of 
SfSj —ni-rij/A is displayed as the thickness of each line, for a 
9x9 system with one hole and J = 0.01. The ferromagnetic 
polaron, with its magnetization direction in the x — y plane, 
is clearly visible in the central region, where the expectation 
value is nearly -0.25. 

We remark that the spin-orientation of the polaron is 
also a physical quantity in the t — J model since the 
Neel order picks out a particular spin orientation (z). 
We might again expect that the polaron spin be oriented 
in the xy-plane, perpendicular to the Neel order vector. 
We may study this issue on a finite lattice by applying a 
small staggered field at the boundary of our lattice. The 
result, shown in Figs. 8 and 9, shows that indeed the po- 
laron spin is oriented in the xy-plane. In Fig. 8, one sees 
that the energy as a function of the total z component 
of spin is a minimum for S z — 0, corresponding to an 
orientation of the polaron in the x — y plane. The back- 
ground antiferromagnetic order is in the z direction. The 
energy cost of rotating the polaron towards the z direc- 
tion is small, however. Once we take S z > 8, the energy 
rises rapidly, indicating that now the polaron is forced to 
change size, rather than simply reorienting. In Fig. 9, we 
show separately the exchange and kinetic energies. The 
interpretation just discussed based on the total energy in 
Fig. 8 is supported strongly by the sharp kinks in the 
energies at S z = 8 in Fig. 9. For S z > 8, the polaron is 
made larger, with a gain in Ek but with a greater cost in 
Em- One sees that for S z < 8, the driving force in orient- 
ing the polaron is the kinetic energy, which is lowest for 
S z = 0. It appears that for ani-y polaron, the hole is 
able to penetrate more readily outside the ferromagnetic 
region, lowering the kinetic energy. 
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FIG. 8. The total energy of a hole, measured relative to 
the undoped system, as a function of J, on an open 9x9 
system, with a staggered magnetic field of ±0.1 applied to all 
edge sites. 
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FIG. 9. (a) For the same system as in Fig. 8, the exchange 
energy of a hole as a function of S z . (b) The kinetic energy 
of the hole. 



We thus conclude that the simple Nagaoka spin po- 
laron picture for a single hole in the t — J model is valid 
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for very small J, less than 0.02 — 0.03, and also is valid for 
the t— J z model. However, analytic estimates for size of 
the polaron are not very accurate in the numerically ac- 
cessible regime. Furthermore, analytic estimates for the 
maximum J at which it a polaron exists are off by about 
a factor of five. For larger J the Brinkman-Rice picture 
of self-retracing excursions into the surrounding antifcr- 
romagnetic region gives energies in fairly good agreement 
with our results. 
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